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In the near infrared and visible bandpass, optical propagation theory conventionally 
assumes that humidity docs not contribute to the effects of atmospheric turbulence on 
optical beams. While this assumption may be reasonable for dry locations, we demonstrate 
in this paper that there is an unequivocal effect due to the presence of humidity upon the 
strength of turbulence parameter, C^, from data collected in the Chesapeake Bay area 
over 100-m length horizontal propagation paths. We describe and apply a novel technique, 
Hilbert Phase Analysis, to the relative humidity, temperature and data to show the 
contribution of the relevant climate variable to as a function of time. © 2008 Optical 
Society of America 

OCIS codes: 010.1300, 010.7060, 030.7060 
1. Introduction 

It has been known for some time^ that the scintillation behaviour of point sources is a measure of the optical 
seeing in the atmosphere. What has been less well understood is the contribution of different environmental 
variables to optical seeing. 

Over the past decade, a great deal of study has been dedicated to clarifying this issue. Comprehensive 
treatments of the theory of wave propagation in random media are given in Tatarskii's seminal works 
Some of the simplest models based on these complex works are well known and available in the literature: 
Huf nagel- Valley P SLC-Day and SLC-NightP These models are commonly used to predict the strength of 
weak clear air turbulence's refractive index structure parameter, C„. 

The underlying assumption is that the index of refraction of air depends solely upon pressure and tem- 
perature at visible and near infrared wavelengths We can verify this from landbound horizontal path 
measurements. These show that during the daylight hours on a clear sunny day, the weak scintillation 
regime trend is dominated by solar insolation, as is illustrated by the left hand plot of Fig. [TJ 

The right hand plot of Fig. [T] illustrates the well known fact that the effects of turbulence do not cease 
after sunset. Here it is less clear as to the predominant contributing factors. The assumption of temperature 
dominance is sensible in the driest of environments, but we will demonstrate from experimental evidence 
that this is an incomplete condition where there is a significant amount of water vapour entering or leaving 
air volume where C„ is measured. 
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1. A. Past work 

A body of work has been reported whose aims have been to compare bulk climate estimates of with 
optically derived readings In those works, the wavelengths of the optical measurements were in the 
mid-infrared (3 to 5 fim), and all showed that humidity along the measurement volume plays a significant 
role in the temporal behaviour of C„. 

In a theoretical study of the structure parameter using bulk measureables over sea ice and snow, Andreas^ 
defined a refractive index scale n* via 

n* = At* + Bq*^ J£- = 1 + (1) 

where i» and q* represent temperature and humidity, scaled such they are constant with height and A and B 
are constants given in Appendix A. The variable n* therefore represents a horizontal refractive index layer 
above the terrestrial surface. He rewrote Eq. ([1]) in terms of the Bowen ratio Bo, the ratio of sensible heat 
flux to latent heat flux, such that 

where c p represents the specific heat of air at constant pressure, L the latent heat of sublimation of ice and 
p is the air density. Since Bo oc this formulation allows the dependence of the refractive index scale 

upon temperature and humidity scales to be easily interpreted. 

We adapt Eq. © to our experimental conditions by setting L to be the latent heat of vapourisation of 
water, and we assume P = 1000 hPa, T = 25 °C and A = 0.93 /an. The result is shown in Fig. (2j from 
this it is clear that when \Bo\ is large, n^f At* is near 1. This means that when sensible heat (heat exchange 
without change of thermodynamic phase) dominates the Bowen ratio, the t* in the denominator of Eq. (TTJ) 
is the principal parameter. 

When \Bo\ is small however, the latent heat flux dominates so q* makes the major contribution to n*. Thus 
there is an interplay between humidity and temperature effects upon the refractive index in the visible/near 
infrared, which vary as a function of the Bowen ratio. The pole at the origin of Fig. [2] is indicative of a 
breakdown of this formulation at very small Bowen ratios. Nevertheless the curve strongly suggests that as 
long as there is significant injection of moisture into the measurement volume of C^, humidity effects become 
non-negligible . 

In this paper we demonstrate that the temporal behaviour of the optically determined turbulence structure 
parameter in the visible and near infrared bandpass at a coastal site is strongly influenced by local 
humidity. Although no Bowen ratio values could be determined from the available measurements, the results 
do support the qualitative interpretation of Andreas' work. 

2. Experiment and Correlogram Analysis 

The C\ and associated weather variable data used in this study were collected over a number of days during 
2003 and 2004 at Chesapeake Beach in Maryland, USA. 

The data were obtained with a commercially available scintillometer (model LOA-004) from Optical 
Scientific Inc, which serves as both a scintillometer and as an optical anemometer for winds transverse to 
the beam paths. The separation between source and receiver was just under 100-m, placing the system in 
the weak scintillation regime. The local weather parameters were determined by a Davis Provantage Plus 
(DP+) weather station. The LOA-004 had a sample rate of 10 seconds, while the DP+ was set at 5 minutes. 

The LOA-004 instrument comprises of a single modulated infrared transmitter whose output is detected 
by two single pixel detectors. The transmitting LED has a bandwidth of 0.92 to 0.96 ^m, while the detector 
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bandwidth is much broader, at 0.65 to 1.0 /^m. The path integrated measurements are determined 
by the LOA instrument by computation from the log-amplitude scintillation (C x (r)) of the two receiving 
signalsP^El] The algorithm for relating C x (r) to C\ is based on an equation for the log-amplitude covariance 
function in Kolmogorov turbulence by Clifford et aZ.P^l 

The C\ data was smoothed with a 60 point (5 minute) rolling average function. We define the morning 
and night portions of a 24 hour period with respect to the measured solar irradiance function, such that we 
exclude the effect of solar insolation from the data in this study. Morning runs from midnight until sunrise, 
while night runs from sunset until 23:59. As reported in Oh et aZFSEED visual inspection of the time series 
data gives the impression that there is an approximate inverse relationship between and relative humidity. 
This can be appreciated in a more quantitative manner by graphing one parameter against the other. 

We chose data sets in which the temperature variations are no more than 16° F and the pressure change 
is at most 15 hPa over the time intervals of interest. The actual range of variation per data set is shown 
in table [T] The data sections were also selected to have no scattering effects due to snow or rain, and the 
wind was northerly (to within approximately ±20°, inflowing from the bay to land). Only a small subset 
eight morning and evening runs, spanning seven days between November 2003 and March 2004, provided 
complete time series in both ambient weather variables and C^. Part of this is shown in Fig. [3J the 
vs. humidity correlograms evidence a negative gradient in all eight runs. Table [1] shows us that the relative 
humidity variation is more strongly reflective of changes in absolute moisture content than of moisture 
holding capacity of the air. 

We conclude that humidity plays a significant part in the behaviour of the refractive index structure 
parameter in the datasets studied. 

3. The Stationarity Problem 

The unsophisticated use of the correlogram tool in Section [2] is a rapid, first order method for inferring the 
statistical influence of one measurable upon another. However it is unsatisfactory as it does not reveal any 
detailed information, such as exactly when the humidity contribution is important with respect to all other 
possible parameters (e.g. temperature) and to what degree its contribution is influential. Cross covariance 
techniques are a possible second order methocP^for extending the analysis, but there are two major stumbling 
blocks. The first is the difficulty in interpreting the meaning of the covariance functions physically and the 
second is the non-stationary nature of the data. In the latter case it is well known that strict stationarity is 
a constraint that is impossible to satisfy on practical grounds, since no detector can cover all possible points 
in phase space. This has motivated us to employ a novel analysis technique to infer the relationship between 
C\ and climate parameters, which we briefly describe in the following section. We leave a more detailed 
study of the limitations of the technique to a separate paper. 

3. A. The Analytic Signal 

GaboiSs] defined the complex analytic signal, namely 

*(£) = X(t) + iY(t) (3) 

— 1 f°° X(t) 
where Y(t) = U[X{t)] = — P / . 1 \ dr 

n J-oo {t - T) 

where represents a Hilbcrt Transform. Through this definition, ^(t) is unique. The Hilbcrt Transform is 
a well known integral transform with a singular kernel (l/(7r(t — t)), t also being a time variable if t is time. 
It is also a Cauchy Principal integral, which we denote by P. Note that the Hilbert Transform preserves the 
norm of the real signal; the difference between Y and X is a ir/2 phase shift for a periodic function. That 
is, the Hilbcrt Transform of cost is sint, and sini is — cost, and the transform of a constant is zero. 
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Eq. ([3]) means that we can write the analytic signal as 



a(t) exp' 



i*(t) 



where a(t) 



y/X*(t)+Y*(t) 



(4) 




which is similar to the well known Fourier expression. We may now determine the instantaneous signal phase, 
Thus we can also calculate the instantaneous frequency u)(t), defined as 



Both and u>(t) can be interpreted as physical measureables, provided certain preconditions are met 
which we describe below. 

3. B. Hilbert Phase Analysis 

There is a problem with the definition of Eq. ([3]); the Hilbert Transform's kernel represents a non-causal 
filter, of infinite support. So if one applies the Hilbert Transform directly to a time varying signal which 
has a non-zero local mean in any subsection, there is a high probability that at least one of a number of 
paradoxes^ will be encountered. This may be best appreciated if we consider a signal in phase space, as in 
the left hand plot of Fig. 2J Here we show the real versus imaginary components of the Hilbert Transform 
of a series of analytic C% measurements. The trajectory of the analytic signal's vector is subject to many 
alterations in the vector start point, norm and phase angle. Any attempt to determine the instantaneous 
frequency is bound to be problematic as we are not able to follow the signal vector's start position over 
time; this starting point obviously does not remain at the origin of the coordinate system. Calculating the 
instantaneous frequency from the Hilbert Transform as is will generate both positive and negative values, 
rendering it physically uninterpretable. 

The paradoxes may be avoided by the application of the Empirical Mode Decomposition (EMD) method 
developed by Huang et aZp^l which we have implemented^^. EMD is a unique and novel method that is 
able to separate an arbitrary real time series into "eigenfunctions" termed Intrinsic Mode Functions or IMFs, 
each of which possesses a structure with well defined instantaneous frequencies, u>(t). The term cigenfunction 
is used suggestively here; we do not mean to imply that the IMFs are eigenfunctions in the strict sense. 

Briefly, the EMD technique consists of (1) the determination of two envelope functions about the time 
series, covering all the local minima and maxima respectively; (2) the computation of the mean of the two 
envelope functions; (3) the removal of the mean of the envelopes from the time series in an iterative manner, 
until the mean is found to be zero; (4) the storing of the resulting mode as an IMF. The IMF is then 
subtracted from the original time series and steps (1) to (4) are repeated, thereby sifting out a family of 
IMFs, stopping only when the resulting mode shows no variation. The final mode represents the overall trend 
of the signal and is not itself an IMF. 

The Hilbert Transform of the IMFs, one of which is shown on the right side of Fig. [4] ensures that 
the analytic signal vector's origin stays fixed and no sudden changes in the direction of u)(t) occur. These 
conditions being satisfied, the instantaneous frequency remains positive and physically meaningful. 

4. Hilbert Phase Analysis 

The Hilbert Phase Analysis (HPA) technique is based on the ideas mentioned in the previous subsections. It 
is clear that a phase angle ($(£)), as well as an amplitude (a(t)), can be found from the Hilbert Transform 
of the IMFs derived from the EMD sifting process. 

The procedure for the HPA data analysis follows a three step protocol: 
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1. Decompose the time series measurements of differing parameters via EMD, obtaining their IMFs and 
trend lines. 

2. Apply the Hilbcrt Transform to the various IMF sets. We discard the trend lines as the reason that 
the original time series fail to be Hilbcrt Transformed in a physically comprehensible manner. 

3. Examine the instantaneous phase angles of the Hilbcrt Transformed IMFs between different parameters 
to infer the dynamics of the physical system. 

4- A. Interpretation of HPA 

In this section we demonstrate that physical effects of a non-linear, non-stationary, time varying system can 
be studied via the sum of all IMF phases, Eja/f ^K*)- 

In Fig. [5] the Ejmf ^(*) * s g ra P ne d with the solar radiation data superimposed for example days. There 
is clearly a gradient change in the phase function at sunrise and sunset. Also notable is the phase jump 
towards a lower gradient whenever the solar radiation function exhibits a drop in amplitude; likewise the 
phase gradient increases with sudden increases in the measured solar radiation. This can be understood as: 
a reduction in energy in the system leads to a lower instantaneous frequency, therefore we see a lower phase 
gradient. A change in energy results in a change in instantaneous frequency, so we sec a modified gradient. 

We conclude that overall physical features of a non-stationary time series can be extracted by inspection 
of the sum of its Hilbcrt Phases of the IMFs. 

5. HPA of C^, humidity and temperature 

5. A. Phase locking between measureables 

To better understand the dependence of <I>c upon humidity and temperature, we consider the difference 
between observable phases, as illustrated in Fig. [7] We define the following difference terms, 

Act = * c ~ $T ' AcH = ^2®c ~^2<f>H (6) 

A CHT = * C 2 

representing the phase differences between and the bulk climate parameter. It is anticipated that the 
controlling climate parameter will be indicated by a near zero difference, A. We reason this because if A 
is constant between C\ and a climate variable then the two datasets must be phase locked. If the mean A 
is zero, then any variations about zero should indicate a synchronization between the datasets. Supporting 
empirical evidence is found upon studying cases wherein the solar insolation influence upon the measured 
is strong, as in Fig. [51 In those plots one can see that the difference between the summed Hilbert Phase 
of C\ and solar radiation (E ~ E ®s) flattens out when the solar function is significant. 

Motivated by this, we determined by linear regression the mean gradients of all the A curves, which we 
list in Table [2j If we assume that the only major contributor should be the local temperature variations, 
Act should be zero or near zero in all cases. This turns out not to have occurred. In fact, Act has no near 
zero value, so temperature variation cannot be said to be phase locked with the refractive index structure 
parameter. 

If humidity has a controlling effect in the visible / near infrared, then it should exhibit a range of values 
near zero for Ach- It only does so for 11/10/03 (am) and 03/28/04 (am). It is noteworthy that 11/09/03 
shows a possible phase lock between the mean of humidity and temperature with C\. This seems to indicate 
that that temperature is vying with humidity for influence over the C\ parameter during that morning. 
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5. B. Further study 

If we define "dominance" to mean that the phase gradient is the value closest to zero of the differences under 
scrutiny then we hnd that humidity is dominant for the 11/10/03, 02/02/04, 03/27/04 and 03/28/04 datasets. 
The 11/03/03 and 04/03/04 (a.m. and p.m.) datasets show that Act < ^ch, indicating the dominant effect 
is temperature, this being extremely strong for the 04/03/04 (pm) dataset and rather weaker for the other 
two. We postulate that the 11/03/03 and 04/03/04 datasets indicate both temperature and humidity are 
contributing to the behaviour of C 2 , with the HPA method possibly showing the proportional contribution 
to C 2 of each climate variable. This is clearly an area to be examined further. 

6. Conclusions 

From experimental data we have shown conclusively that humidity plays a significant part in the visible/near 
infrared measure of C 2 in a coastal environment. This is in qualitative agreement with Andreas' model and 
is a natural extension of the results from the mid infrared. 

Furthermore we have explored a new technique, Hilbcrt Phase Analysis, with which to study this physical 
phenomenon. In overall terms, the HPA method is in agreement with the correlogram results. Phase locking 
between data is an unexpected result that needs further examination. We have found that the HPA technique 
described here is very promising and is likely to provide much more information about the changes to a 
physical system than traditional methods. 

Appendix A: Definitions 

The constants A and B of Eq. Q]in the 0.36 to 3 (im wavelength region are 

A = -10" 6 mi(A)(P/T 2 ) (Al) 
B = 4.6150 x -10~ 6 [to 2 (A) -mi (A)]. 

P and T are the ensemble average air pressure and temperature, respectively. The functions mi and mi are 
polynomials in wavenumber given by 

„o W ,o. 6839.397 45.473 

mi (A) = 23.7134 H , , + , , xn A2 

u ; 130 -(1/A) 2 38.9 -(1/A) 2 v ; 

ma (A) = 64.8731 + 0.58058(1/A) 2 - 0.0071150(1/A) 4 + 0.0008851(1/A) 6 
with A's units in micrometres. 
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List of Figure Captions 

Fig. 1. The effects of solar insolation on C 2 as measured on March 9, 2006 in Puerto Rico. Solar radiation 
measurements are superimposed at an arbitrary scale on top of the C 2 data. (Left) The full 24 hour period. 
The vertical axis values are of order 10 -12 to -2 / 3 . (Right) Close up of the evening period. The vertical axis 
values are of order 10~ 14 uT 2 / 3 . The C* 2 data were obtained at an urban site 1.75-km from the sea. These 
values are in agreement with measurements of seeing over the sea by alternative meansP^ 
Fig. 2. Graph of n^/(At^) vs Bowen ratio for A = 0.93 /im. 

Fig. 3. Example correlograms of C 2 and Relative Humidity in the absence of solar insolation. The upper and 
lower bounds indicate the 50% confidence level. The C 2 magnitudes are 10 -15 m -2 / 3 , in agreement with 
measurements of seeing over the sea by alternative means 

Fig. 4. (left) Hilbert Phase Space plot of the trajectory of the C 2 signal vector of Fig. [T] The signal vector's 

start point drifts arbitrarily around the phase space, making the determination of a physically reasonable 

instantaneous frequency impossible over the whole path, (right) Hilbert Phase Space plot of the trajectory 

of a single IMF derived from the C 2 signal. The IMF vector's start point is stable and its trajectory does 

not change direction, so a positive instantaneous frequency can be determined at all points. 

Fig. 5. Sum of all Hilbert Phases of the measured Solar Radiation IMFs Q2 $5) for the following days: Feb 

2, Mar 27, Mar 28 and Apr 3, 2004. The Solar Radiation is superimposed at an arbitrary scale. 

Fig. 6. Plots of Hilbert Phase differences between C 2 and solar radiation for 03/27/04 and 03/28/04. The 

solar radiation function is superimposed at an arbitrary scale on each graph. Note the flattening out of the 

phase difference function during the daylight hours. 

Fig. 7. Difference plots for the morning of November 10, 2003. The dotted line is a linear regression, estimating 
the mean phase gradient. The top left graph shows a phase lock between C 2 and relative humidity. 
Table 1. Mean and range of bulk parameters measurements and the maximum range of specific humidity 
(AqJ. 

Table 2. Linear regression line gradients of the Phase Differences. 
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Fig. 1. The effects of solar insolation on C 2 as measured on March 9, 2006 in Puerto Rico. 
Solar radiation measurements are superimposed at an arbitrary scale on top of the C 2 data. 
(Left) The full 24 hour period. The vertical axis values are of order 10 -12 m -2 ' 3 . (Right) 
Close up of the evening period. The vertical axis values are of order 10~ 14 m~ 2 / 3 . The C 2 
data were obtained at an urban site 1.75-km from the sea. 
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Fig. 2. Graph of n*/ (At*) vs Bowen ration for A = 0.93 ^m. 
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Tabic 1. Mean and range of bulk parameters measurements and the maximum range of 
specific humidity (Aq s ). 



Date 


RH 

(%) 


T 

(°F) 


P 
(hPa) 




A T 


A P 


A Qs 
(g/kg) 


11/03/03 (pm) 


83.6 


65.2 


1015.3 


15 


2 


0.4 


2.2 


11/09/03 (am) 


49.3 


38.8 


1030.9 


11 


3.1 


1.2 


0.5 


11/10/03 (am) 


76.4 


46.8 


1027.6 


11 


8 


1.8 


0.6 


02/02/04 (am) 


91.4 


29.9 


1022.9 


16 


2.5 


6.6 


0.8 


03/27/04 (pm) 


76.4 


46.8 


1027.6 


11 


8 


1.8 


1.7 


03/28/04 (am) 


71 


44.5 


1027.6 


16 


0.7 


0.6 


0.9 


04/03/04 (am) 


76.1 


46.6 


997.4 


18 


5.9 


4.9 


0.4 


04/03/04 (pm) 


43 


40.1 


1004.5 


14 


5.6 


2.1 


1.1 
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Fig. 3. Example correlograms of and Relative Humidity in the absence of solar insolation. 
The upper and lower bounds indicate the 50% confidence level. The C\ magnitudes are 10~ 15 
to -2 / 3 , in agreement with measurements of seeing over the sea by alternative meansP^l 
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Fig. 4. (left) Hilbcrt Phase Space plot of the trajectory of the signal vector of figure 
[TJ The signal vector's start point drifts arbitrarily around the phase space, making the 
determination of a physically reasonable instantaneous frequency impossible over the whole 
path, (right) Hilbert Phase Space plot of the trajectory of a single IMF derived from the 
signal. The IMF vector's start point is stable and its trajectory does not change direction, 
so a positive instantaneous frequency can be determined at all points. 
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Fig. 5. Sum of all Hilbert Phases of the measured Solar Radiation IMFs 3>s) for the 
following days: Feb 2, Mar 27, Mar 28 and Apr 3, 2004. The Solar Radiation is superimposed 
at an arbitrary scale. 
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Z * - Z 4> : 03/28/2004 
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24-hour time 24-hour time 

Fig. 6. Plots of Hilbert Phase differences between and solar radiation for 03/27/04 and 
03/28/04. The solar radiation function is superimposed at an arbitrary scale on each graph. 
Note the flattening out of the phase difference function during the daylight hours. 



Table 2. Linear regression line gradients of the Phase Differences. 



Date 




Act 


^CHT 


11/03/03 (pm) 


-9.7(4) 


-8.0(0) 


-8.8(7) 


11/09/03 (am) 


-11.2(4) 


9.7(7) 


-0.7(4) 


11/10/03 (am) 


-0.4(1) 


12.7(9) 


6.1(9) 


02/02/04 (am) 


8.5(1) 


-12.3(0) 


-1.8(9) 


03/27/04 (pm) 


3.7(6) 


28.1(7) 


15.9(7) 


03/28/04 (am) 


0.0(2) 


-3.3(1) 


-1.6(5) 


04/03/04 (am) 


6.8(1) 


6.6(8) 


6.7(4) 


04/03/04 (pm) 


-27.3(3) 


-2.3(3) 


-14.8(3) 
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Fig. 7. Difference plots for the morning of November 10, 2003. The dotted line is a linear 
regression, estimating the mean phase gradient. The top left graph shows a phase lock 
between C\ and relative humidity. 
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